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Abstract. Among their competitors, projection deptli and its induced esti- 
mators are very favorable because they can enjoy very high breakdown point 
robustness without having to pay the price of low efficiency, meanwhile provid- 
ing a promising center-outward ordering of multi-dimensional data. However, 
their further applications have been severely hindered due to their computa- 
tional challenge in practice. In this paper, we derive a simple form of the pro- 
jection depth function, when (/i,cr) = (Med, MAD). This simple form enables 
us to extend t he existing result of point-wise exact computation of projection 
depth (PD) of Zuo and Lail ( 201ll ) to depth contours and median for bivariate 
data. 
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1 Introduction 

To generalize order-related univariate statistical methods, depth functions have emerged as 
powerful tools for nonparametric multivariate analysis with the ability to provide a center- 
outward ordering of the multivariate observations. Points deep inside a data cloud get higher 
depth and those on the outskirts receive lower depth. Such depth induced ordering enables one 
to develop favorable new robust e stimators of multivariate location and scatter matrix. Since 
Tukey's introduction ( Tukey . I975I ). depth functions have gained much attention in the last two 
decades . Num erous depth notat ions have been introduced. T o name a few, halfspace dept h 
( Tukey . I975I ). si mplicial depth (Liu . I99C ). regres s ion depth ( Rousseeuw and Hubert . I999I ). 
projection depth (|Lh] . ll992l : rZuo and Serflin3 . l200d : Iz^ . booi ). 

Zuo and Serfling (|200d l and S (|2003l l found that among all the examined depth notions, 
projection depth is one of the favorite, enjoying very desirable properties. Furthermore, pro- 
jection depth induced robust estimators, such as projection depth weighted means and median, 
can possess a very high breakdown point as well as high relative efficiency with appropriate 
choices of univ ariate loca. t ion and scale est imators, serving as very favorable alternatives to the 
regular mean ( Zuol . 2003 : Zuo et al . 20041 ). In fact, projection depth weighted rneans i i iclude 
as a special cas e the famous Donoho-S tahel estimator ( Stahel . 1981; Donohol . 1982 : Tyler . 1994; 
Maronna et al . 19951 : Zuo et al. . 20041 ). the latter is the first constructed location estimator in 
high dimensions enjoying high breakdown point robustness and affine equivariance, while the 
projection depth median has the highe st breakdow n point among all the existing affine equiv- 
ariant multivariate location estimators ( Zuol . 20031 ) . 



'Corresponding author's email: zuo@msu.edu, tel: 001-517-432-5413. 
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However, the further prevalence of projection depth and its induced estimators is severely 
hindered by their computational intensity and intimidation. The computation of projection 
depth seems intractable since it involves supremum over infinitely many direction vect ors. There 
were only approximating algorithms in the last three decades until Zuo and Lai ( 201ll ). in which 
they proved that there is no need to calculate the supremum over infinitely many direction 
vectors in the bivariate data when the outlyingness function uses the very popular choice (Med, 
MAD) as the univariate location and scale pair. An exact algorithm for projection depth and 
its weighted mean, i.e. Donoho-Stahel estimator, was also constructed in that paper. 

In the current paper, we further generalize the ir idea to the higher dimensional cases by 
utilizing linear fractional functionals programming (' SwarupI . EoG^ '). That is, we find that, with 
the choice of (Med, MAD), we only need to calculate the supremum over a finite number of 
direction vectors for p > 2. Furthermore, these direction vectors are x-free and depend only 
on the data cloud. Therefore, we derive a simple form of the projection depth function, and 
are able to compute the bivariate projection depth co ntours and median very conveniently 
through linear programming based on the procedure of Zuo and Lai ( 201ll ). It is found that 
sample projection depth contours are polyhedral under some mild conditions. Furthermore, 
it is noteworthy that the computational methods discussed in this paper have no limitation 
on the dimension and therefore could possible be im plemented to spaces with p > 2, as 
well as to the modified projection depth ( Siman . 201ll ) in a more general multidimensional 



regression context. The corresponding programs are available from the authors (zuo@msu.edu 
or csuliuxh912@gmail.com). 

The rest of the paper is organized as follows. Section [2] provides the definitions of the 
projection depth contour and projection median. Section [3] presents the main idea of how to 
get a simple form of the projection depth function. Section [J] discusses the exact computational 
issue of the projection depth contour and projection median by linear programming. While some 
numerical examples are given in Section [5l 



2 Definition 



For a given distribution F on M^, let /i(-F) be translation equivariant and scale invariant, 
and (j{F) be translation invariant and scale equivariant. Define the outlyingness of a point 
X € B P (p > 1) with respect to the distribution F of the random variable X G BP as (see (jZud . 
20031 ) and references therein) 

0{x,F) = sup |Q(u,x,F)| (1) 



where Q(n, X, F) = {u'^x — fi{Fu))/cr{Fu), if u'^x — fi{Fu) = cr{Fu) = 0, then define Q{u, x, F) 
Fu is the distribution of u'^X, which is the projection of X onto the unit vector u. 



0. 



Throughout this paper, we select the very popular robust choice of fx and a: the median (Med) 
and the median absolute deviation (MAD). Based on definition ([H), the projection depth of any 
given poiri t x wit h respective to F, PD{x, F), can then be defined as (|Lh] . (l992l :l IZuo and Serfiina . 
200d : S, l2003l ^ 

PD{x,F) = l/(l + 0(x,F)) 

With the outlyingness function and projection depth function defined ab ove, we then define 
the projection depth median (PM) and contours (PC) as follows ( Zuo . 20031 ) 



arg max PD(x, F), 

xeRp 



PM{F) 

PC{a, F) = {xeRP : PD{x, F) 



a} 
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where < a < a* = sup^jg^p PD{x, F). 

For a given sample A*" = {Xi,X2, ■ ■ ■ ,Xn} from X, let Fn be the corresponding empirical 
distribution. By simply replacing F by F„ in PM{F) and PC {a, F), we can obtain their sample 
version: PM{Fn) and PC(a, Fn). Without confusion, we use X"" and Fn interchangeably in what 
follows. Furthermore, by noting the fact that for the choice of (Med, MAD), Q(u,x, X"') in ([T|) 
is odd with respect to u, we drop the absolute value sign existing in definition ([1]), and consider 

0(3;, A"") = sup 

||n||=l 

instead in what follows, where 

u'^x - Med(u^A'") 



Q{u,x,X'') 



MABiu'^X'^) 



where u'^x denotes the projection of x onto the unit vector u, and u'^X^ = {u^ Xi,u^ X2, • • • , u^Xn}- 
Let < Z(2) <, • • • , < ^(n) be the order statistics based on the univariate random variables 
Z" = {Zi,Z2,--- then 

Med(Z") = ^(L("+i)/2J) +^(L(n+2)/2j) 

MAD(Z") = Med{|Zi - Med(Z")|, i = 1, 2, • • • , n}, 
where [-J is the floor function. 



3 The main idea 



Note that, for any given sample A'", the tasks of computing both PM(X^) and PC{a, X"^) 
mainly involve O^XjX""), i.e. 

PM^X"") = arg min 0(x, Af"), 

x^RP 

PC{a, Af") = {xeRP : 0{x, X"") = /?} , 

where (3 = \/a — 1. Thus, let's first focus on the computation of 0{x,X^). Without loss of 
generality, in what follows, we assume X ^ to be in general position , which is commonly used in 
most existing literature; see for example iDonoho and Gaskol (jiggd l. 

By the idea of a circular sequence dEdelsbrunnei] . [lOS^) (see also IPvckerhofj jiooS); ICascosI 
( 20071 )). for any given unit vector v € S = {u £ BP : \\u\\ = 1}, there must exist two permuta- 
tions, say (n,i2,-- - ,in) and (ji,j2,-- - ,Jn), of (1,2,- •• ,n) such that 



v-'Xi, <v^Xi, <---<v^Xi^, 



Y < Y < ■ ■ ■ < Y 



v'^Xj^ — Med{v'^X'^)\, {1 < I < n). There is a small non-empty set J\f{v) C S of v 
such that these hold true for any v G J\f{v). This implies that the whole unit sphere <S can be 



where Yj^ 
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covered completely by at most N (with order )) non-empty fragments 5^ = {f € 5: v 

satisfy constraint conditions Qk} with Qk being 

Y —Y > n 



jk.n jk.n — 1 — ' 

I ' 11^11 = 1, 

for some fixed permutations {ik,i,ik,2, ■ ■ ■ ,ik,n) and (ifc,i,Jfc,2, • • • ,jk,n) of (1, • • • where 1 < 
k<N. 

Note that different fragments Sk {k = 1,2, ■■■ , N) are connected and overlapped with each 
other only on the boundaries. Thus, to calculate 0{x,X^), it is sufficient to calculate 

0(x,;f") = max Ofc(x,;f") 

1<A;<A'' 

with 

Ok{x, X"^) = sup Q(n, X, ;f"). (2) 



Furthermore, from the definition and property of Sk, it is easy to see that, for any u ^ Sk, 
the outlyingness function Q{u, x, X^) can be simplified to 

u'^ix — Xi, ) 

I Jk,m ''k,m I 

if n is odd with m = (n + l)/2, otherwise 

Q(u,.,r') = -. ^"'(--x-^^;) ^, (4, 



\U' 



with m = n/2, m* = m + 1 and Xk^a = (^i^,™ + ^i^^^, )/2. 

Remark 1. Based on the assumption that X"^ are in general position, the denominators in 
the above two formulas will not be for all u (z Sk, since they are actually equal to M.ADju'^ X"'), 
and g reater than under such an assumption when n > 2p; see the proof of Theorem 3.4 in lZuo 
(I2OO3I ). 

By ([3]) and ([!]), we obtain the following proposition. 

Proposition 1. Assume X"' is in general position. Then for any given A; (1 < A; < N), the 
optimization problem ([2]) is equivalent to 

0,(x,;f-) = sup^, (5) 

z ^k^ 

subject to 

Akz > (6) 
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where Ck, dt and will be specified in the Appendix. Here b > means that b is component- 
wise non- negative if b is a vector, i.e. for any component bj, we have bj > 0. 



([5]) with constraint cond itions (E l) is ty pically a l inear fraction al functionals programming 
problem. By theorem 1 of ISwarupl ( 19621 ) (see also ISima 3 (1201 ih (p. 950) for more general 



discussion), it is easy to show that the maximum of c^z/d^z will only occur at the basic feasible 
solution of ([6]). Note that the number of fragments is limited (at most N). Thus, we have 

Theorem 1. Suppose that the choice of location and scale measures of projection depth 
function is the pair (Med, MAD). Then the number of direction vectors needed to compute the 
projection depth exactly is finite. Furthermore, these direction vectors only depend on the data 
cloud A'". 

Remark 2. The idea of dividii ig the unit sphere S into fragments by applying Med and 



MAD sequenc es was first used in IZuo and Lail (j201ll ) for computing the bivariate projection 
depth; see al so iPaindavein and Siman (|2011bH d^ for other similar applications. Here we extenci 
the result of IZuo and Lail (|201ll ) to BP {p > 2). That is, one could compute PD in exactly 
by only considering a finite number of direction vectors. Furthermore, the x-free property of 
these direction vectors can bring convenience to the computation of PD(x, X"") for any x, since 
we only need to search the direction vectors once. 

4 Exact computation of PM(A:'") and PC{a, A'") 

From the discussion above, we can obtain the two following observations, namely, for any 
given X, 

• the way to divide sphere S into fragments Sk {k = 1,2, ■ ■ ■ , N) is fixed, i.e. j;-free, as long 
as the data cloud X"" is fixed. 

• there is no need to calculate Q{u, x, Af") over an infinite number of direction vectors. It is 
enough to calculate it for, say, ui,U2, - • • ,um- 

Based on the discussion and two observations above, we therefore can re-express the outly- 
ingness function 0{x, X"^) as follows 



0{x,X'') = max 



u\x - Med{ulX'') u\x - Med(u5;f") u\iX - MaAiul^X' 



\^kV>{u\X-^) ' MAD(u5A'") ' ' yVKD{u\^X'^) 



where \ui\^L^^ are some p-dimensional vectors depending only on the data cloud <Y". For the sake 
of convenience, hereafter we write gi{x) = ajx — bi {i = 1,2, ■■ ■ ,M), where = ]y[AD(n^A'") ^^ 

and h- - ^^'^(^^^^ 



MAD{m[A'") • 

Obviously, 0{x,X'^) = maxi<i<M{gi{x)} is in fact a piece-wise linear convex function with 
respect to x for the given data cloud A^". Therefore, its minimizers can be found by using 
common linear programming methods by solving the problem 

s = min t 

z 

subject to 

t>gi{x), i = l,2,--- ,M, 

where z = {t, x'^Y . This kind of problem can be solved by some common solver such as linprog.m 
in Matlab. Let zq = {tQ,XQy be a final solution of this problem. Then, it is easy to show that 
xq is one of the deepest points with depth value PD{xq, Af") = 1/(1 -|- to) = «*• 
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Given the nature of the maximum piece- wise hnear convex function maxi<j<jvf gi{x), there is 
either a single minimizer or a convex polyhedral set of minimizers. Then there naturally comes 
a question, namely, after obtaining the value a* , how to get all of these vertexes? Note that the 
projection median is a specific case of the projection depth contour. Therefore, let's focus now 
on the computation of projection depth contours. For any gi ven < a < a* , the projection 
depth contour is the boundary of the projection depth region 



{x G RP 
{x € RP 
{x G RP 



PDix^X"^) > a} 

0(x,^")</3} 

gi{x)< f3,i = 1,2,- ■■ ,M}, 



Typically, the regions constrained by linear inequalities such as 

g^{x)<(3,^ = l,2,■■■ ,M 



(7) 



are polytopes. Therefore, the bound ary of PR(a, X^), i .e. PC {a, -Y"), could be easily found by 
employing procedures such as q hull ( Barber et al. j_ 19961 ) based on the dual relationship betwee n 
vertex and facet enumeration (jBremner et aU . 1 19981 ) ; see also iPaindavein and Siman (|2nilal V 
In Matlab, these kinds of tasks can be fulfilled by the function con2vert.m, which was developed 
by Michael Kleder, and can be downloaded from Matlab Central File Exchange. 

However, in many practical applications, the number M (with order 0((p)^)) of the direction 
vectors may be very large. When M is too large, it is difficult to obtain the boundary of the region 
formed by ([7]) by using some of the aforementioned procedures such as con2vert.m, since they 
involve solving some large generalized inverse matrices. Therefore, it is important to eliminate 
some redundant constraints before computing PC{a, X^) for too large M. 

Note that, for any given a (0 < a < a*), the number of the non-redundant constraints in 
([7]) is much small compared to M, which implies that numerous inequalities in ([7]) could be 
eliminated during the computation of the a-contour. In fact, it is not difficult to show that 

Pi2(a,;f") = CinC2n---nc„ 

where = {x : gj{x) < = i^, -|- 1, • • • ,ik+i} with 1 < A; < s — 1 and 1 = ii < 12 < 
■ ■ ■ < is = M, and that Ck is non-empty if mini<Ks_,i{f/_|_i — i/} ^ p. Then, when M is large, 
a procedure for exactly computing PC{a, X"^) is: (1), to find the non-redundant constraints 
in each at first, then (2), to use all of these non-redundant constrains together to compute 
PC{a,X'^). 

With the vertexes in hand, some common graphical packages could be employed to visualize 
these contours very easily in spaces of p = 2, 3. It is noteworthy that, although all the methods 
discussed above could possible be implemented to spaces with p>2 theoretically, feasible exact 
algor ithm for computing the projection depth exists now only in the bivariate data (IZuo and Lail . 
201 ll ). Therefore, we can only provide some exactly results about the bivariate projection depth 
contours and projection median i n the current paper . All the direction vectors {ui, • • • , um} are 
found by using the procedure of Zuo and Lai ( 201ll ). The x-free property of these vectors can 
be proved similarly by using the linear fractional functionals programming as Proposition 1. 



5 Numerical analysis 

In order to gain more insight into the sample version of projection depth contours and pro- 
jection median, we construct some numerical examples in this section. 
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5.1 Simulation results 



To illustrate the robustness and the shape of the bivariate projection depth contours and 
median, we present two examples as follows. The data are mainly generated from the normal 
distribution, but contain a few outliers. 

Example 1. We first generate 60 samples X = {Xi,X2y from the normal distribution 
A^(0,/2), and then disturb these samples randomly by replacing their first components by 6 
with probability 0.05. 

Example 2. We first generate 400 samples X = {Xi,X2Y from the normal distribution 
A^(0, So), and then disturb these samples randomly by replacing their first components by 6 
with probability 0.10, where 

'1 0.5^ 
0.5 1 



So 



For the sake of comparison, the population versions of PC{a,X) corresponding to these two 
examples are also provided here, and plotted according to the formula 



V' 



X2j o? 



which was developed in lZuo (I2OO3I I. Here C^v = ^>-i(f ) « 0.6744898, S denotes the covariance 



matrix of normally distribution X, namely, I2 in example 1 and Sq in example 2. The population 
versions PC{a, X) of example 1 and 2 are given in Figured] and [3l respectively, while the sample 
versions PC{a, X^) are plotted in Figure [2] and SJ respectively. 

Comparing the figures of sample versions with those of population, it is ready to see that 

• The sample versions PC{a, X"') are roughly elliptical, similar to the shapes of the pop- 
ulation versions PC{a, X). Furthermore, it is remarkable that they are resistant to the 
notorious vertical outliers in Figure [2] and HI These res ults confirm that PC{a, X"^) are 
very robust, and could capture the feature of PC{a,X) (Zuo, 20031 ). even when there are 
a few outliers. 

• It is very interesting to note that the shapes of the projection depth contours are also poly- 
lateral due to t he properties of the fu r iction maxi<j<j\,/ (7i(x), similar to that of halfspace 



depth contours (jRuts and RousseeuwI . 119961 ). On the other hand, unlike halfspace depth 



contours, PC{a, X'^) does not need to pass through the observations. 
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X 

Figure 1: Population version PC{a,X) for Example 1 with a = 0.1,0.2, • • • ,0.9 (from out to 
center) . 
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Figure 2: Sample version PC{a, X'^) for a = 0.1,0.2,- •■ ,0.7 (from out to center). Here the 
most inner big point denotes the projection median. 
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Figure 3: Population version PC{a,X) for Example 2 with a = 0.1,0.2, • • • ,0.9 (from out to 
center) 

8 I , , 1 1 , , , 1 



> - 




_8l , , , , , , , 1 

-8 -6 -4 -2 2 4 6 8 

X 

Figure 4: Sample version PC{a, X"') for a = 0.1,0.2, • • • ,0.8 (from out to center). Similar to 
Figure [21 the most inner big point is also the computed projection median. 



Furthermore, to gain more information about the shape of the projection depth contours, we 
also provide some other examples in the following. The sample size we used is 2500. Here Figure [5] 
reports the projection depth contours PC {a, <Y") corresponding to the uniform distribution over 
the triangle with its vertexes being (0, 0), (0, 1) and (1, 1). Figure [6] corresponds to the uniform 
distribution over region [0,1] x [0,1]. Figure [7] gives the contours of the normal distribution 
N{0,l2)- While Figure [8] provides the contours of 5X + {1 - 6)Y, where X ~ iV(/xi,/2), Y ~ 
-^(^21-^2)1 fJ-i = (—2,-2), fi2 = (2,2), and 5 is a discrete random variable with probability 0.5 
taking value 1 and 0.5 taking value 0. Here X, Y and d are all independently disturbed. From 
these figures, we can see that the projection contours are also polylateral and convex. 
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■-1 -0.5 0.5 1 1.5 2 2.5 
X 



Figure 5: Projection depth contours of the uniform distribution over the triangle formed by 
vertexes: (0,0), (0, 1) and (1, 1), where a = 0.1,0.2, ■ ■ ■ ,0.7 (from out to center). 
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Figure 6: Projection depth contours of the uniform distribution over [0, 1] x [0, 1], where a = 
0.1, 0.2, • • • , 0.9 (from out to center). 
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Figure 7: Projection depth contours of the bivariate standard normal distribution A^(0,/2), 
where a = 0.1, 0.2, • • • , 0.9 (from out to center). 
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Figure 8: Projection depth contours of the random vector bX + (1 — (5)y, where a = 
0.2, 0.3, • • • , 0.9 (from out to center). 



5.2 Real data example 



Here a real data example is presented to illustrate the performance of projection depth con- 
tours. 

Tabled] here is taken from Table 7 of Rousseeuw and Lerov ( 19871 ) (p. 57). Total 28 animals' 
brain weight (in grams) and body weight (in kilograms) are presented in this table. Before the 
analysi s, logarithmic transformation was taken for the sake of convenience. According to the re- 
sults of lRousseeuw and LeroyI (Il987l ). there are five cases considered as outlying, i.e. diplodocus, 
human, triceratops, rhesus monkey and brachiosaurus. Among them, the most severe cases are 
diplodocus, triceratops and brachiosaurus. In fact, these three cases are referred to as dinosaurs 
because they possess a small brain as compared with a heavy body (see Table [T]) and their highly 
negative residuals can lead to a low slope for the least squares fit. For the remaining two cases, 
although their actual brain weights are higher than those predicted by the linear model, they 
are not worse than the three previous cases since they do not obey the same trend as that one 
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Table 1: Body and Brain Weight for 28 Animals ( Rousseeuw and Lerov . 19871 ) . 



Index 
i 



Species 



Body Weight 



Brain Weight 
Yi 



1 


Mountain beaver 


1.350 


8.100 


2 


Cow 


465.000 


423.000 


3 


Gray wolf 


36.330 


119.500 


4 


Goat 


27.660 


115.000 


5 


Guinea pig 


1.040 


5.500 


6 


Diplodocus 


11700.000 


50.000 


7 


Asian elephant 


2547.000 


4603.000 


8 


Donkey 


187.100 


419.000 


9 


Horse 


521.000 


655.000 


10 


Potar monkey 


10.000 


115.000 


11 


Cat 


3.300 


25.600 


12 


Giraffe 


529.000 


680.000 


13 


Gorilla 


207.000 


406.000 


14 


Human 


62.000 


1320.000 


15 


African elephant 


6654.000 


5712.000 


16 


Triceratops 


9400.000 


70.000 


17 


Rhesus monkey 


6.800 


179.000 


18 


Kangaroo 


35.000 


56.000 


19 


Hamster 


0.120 


1.000 


20 


Mouse 


0.023 


0.400 


21 


Rabbit 


2.500 


12.100 


22 


Sheep 


55.500 


175.000 


23 


Jaguar 


100.000 


157.000 


24 


Chimpanzee 


52.160 


440.000 


25 


Brachiosaurus 


87000.000 


154.500 


26 


Rat 


0.280 


1.900 


27 


Mole 


0.122 


3.000 


28 


Pig 


192.000 


180.000 



followed by the majority of the data. 

We plot the projection depth contours in Figure El where the green point is the projection 
median with depth value 0.73257, five labeled points denote the outliers mentioned above with 
the points 1-3 corresponding to the case of diplodocus, triceratops and brachiosaurus and 4-5 
corresponding to those of human and rhesus monkey. 8 contours are plotted there. From Figure 
[9l we can see that all of these three dinosaurs lie outside the contour 0.1, while points 4-5 lie be- 



tween the contours 0.1 and 0.15. These results are consistent with those of lRousseeuw and Lerov 



(1983), implying that projection depth contours can capture the structures of the objective data 



and identify outliers. 

Furthermore, it is worth mentioning that the shape of these plotted contours is not af- 
fected by a few atypical points at the outskirts of the cloud, namely, both the inner and outer 
depth contours are roughly elliptical, unlike those of halfspace depth contours (see Figure \W\i 
( Ruts and Rousseeuw . 19961 ). This is the most outstanding difference between projection and 



halfspace dep th contour s, confirming the more robustness property of projection depth and its 
contours (see Zuol f 20041 )). 
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Figure 9: Projection depth contours with a = 0.7, 0.6, • • • , 0.2, 0.15, 0.1 (from center to out). 
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logorithmic body weight 

Figure 10: Halfspace depth contours ( Ruts and Rousseeuw . 19961 ). 
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Appendix Proofs of main results 



Proof of Proposition 1. Here, without loss of generality, we prove only the odd n case. 
Note that, for any u S Qfc, we have 



u^Xi^ < u^Xi^ < ■■■< u^X,^ <■■■< u^Xi^ 



according to the definition of Qk- This implies that 



\u^{X, - Xi^ , 



-u^ {Xi - Xi^^^^), if i e {ik,i,ik,2, - ■ ■ ,ik,m-i}, 

U^{Xi — Xii, if i ^ {ik,m,ik,m+l, ■ ■ ■ ,ik,n}- 



That is, we can remove the absolute value signs of Yj^ based on the order information existing 
in the permutation {ik,i,ik,2, • • • ,ik,n)- Therefore, for any u E Qk, ^ can be further simplified 
to 

Q(u,x,;f" 



where Ck = x - Xi^^rn and dk = Sk{jk,m) ■ (Xj^^m - Xi^^rn) (f <k < N), with 

— I if ^ £ {ik,l,ik,2, ■ ■ ■ , ik,m-l}, 

Next, note that, for any positive A and z = Ait, it holds that cj^z/dj^z = cj^u/d'j^u and bj^z > 
if bj,u > 0. Then, ([3]) and constrain condition Qk lead to 



Ok{x,X"-) = sup 



dlw 



(8) 



subject to 



Akz > 

1 <k < N, with Ak = {^ll), where 

Aki = {XJ^^^ - XI^^^,XJ^^^ - XI 

and 



X'^ X'^ 1^ 



I ^^ijK2)-{Xl.^-XlJ-Sk{jk,i)-{Xl,-XlJ \ 
SkUk,3) ■ {XJ^, - X^J - SkUk,2) ■ {X^, - X^J 



^k2 



\skUk,n) ■ {Xl„ - XIJ - Sk{3k,n-l) ■ (XJ^,„_, - XI J) 

This completes the proof of Proposition 1. 
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